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Abstract 

Finding an optimal alignment connecting two end-points in a speci¬ 
fied corridor is a complex problem that requires solving three interrelated 
sub-problems, namely the horizontal alignment, vertical alignment and 
earthwork optimization problems. In this research, we developed a novel 
bi-level optimization model combining those three problems. In the outer 
level of the model, we optimize the horizontal alignment and in the inner 
level of the model a vertical alignment optimization problem considering 
earthwork allocation is solved for a fixed horizontal alignment. Derivative- 
free optimization algorithms are used to solve the outer problem. The 
result of our model gives an optimal horizontal alignment in the form of 
a linear-circular curve and an optimal vertical alignment in the form of 
a quadratic spline. Our model is tested on real-life data. The numeri¬ 
cal results show that our approach improves the road alignment designed 
by civil engineers by 27% on average, resulting in potentially millions of 
dollars of savings. 


1 Introduction 

Road design optimization is the problem of finding a curve that minimizes the 
construction cost while satisfying all of the desired design specifications. The 
problem is usually divided into three interrelated sub-problems [1]: horizontal 
alignment optimization, vertical alignment optimization, and earthwork opti¬ 
mization. In typical road design, first a horizontal alignment is proposed and 
then an associated vertical alignment is optimized considering the earthwork 
allocations and vertical alignment design constraints. 

For a fixed horizontal alignment, the optimization of vertical alignment is 
a well studied problem [161 0 uni 0 US]. Many of the past approaches have 
modeled the optimization of vertical alignment by way of mixed integer linear 
programming. This creates a complicated, but deterministic optimization prob¬ 
lem that is generally solvable using modern MILP solvers (assuming reasonable 
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road lengths and time allowances) [TS]. This implies that given a proposed hori¬ 
zontal alignment, it is possible to evaluate the quality of that alignment in terms 
of the optimal cost vertical alignment. Which further suggests that it should be 
possible to have a computer evaluate and seek an optimal horizontal alignment 
(in terms of the minimal vertical alignment construction cost). In this paper, 
we demonstrate the practicality of this idea, and further demonstrate the value 
of this approach in terms of cost savings for the final road design. 

To do this, we formulate the horizontal alignment optimization problem as 
a bi-level optimization problem. In the inner level, for a fixed horizontal, a 
vertical alignment optimization problem is solved using the mixed integer linear 
programming (MILP) from [T3] (see Appendix [A| , which builds on min [13 
El] and provides a global optimum. The outer level of the problem is solved 
using a derivative-free optimization algorithm, which gives a local optimum, 
with the starting alignment being the best one produced by a civil engineer. 

The article is organized as follows: Subsection 1.1 overviews some of the past 
research in road design. Section 2 describes basic terminology, Section 3 explains 
the geometric specifications of a horizontal alignment. Section 4 describes our 
proposed horizontal alignment optimization models in detail. Section 5 describes 
the derivative-free optimization solvers we used. Section 6 reports the numerical 
results for the test problems, and Section 7 summarizes the contributions and 
highlights some future works. 

1.1 Past Research in Road Design 

As mentioned, road design is commonly divided into three interrelated sub¬ 
problems: earthwork optimization, vertical alignment optimization, and hor¬ 
izontal alignment optimization. Each problem relies on the solution to the 
sub-problem proceeding it (i.e., vertical alignment optimization requires solving 
earthwork optimization, and horizontal alignment optimization requires solving 
vertical alignment optimization). As such, there is a clear hierarchy in terms of 
problem difficulty. 

Earthwork optimization is perhaps the most established of the three sub¬ 
problems. Many studies investigated earthwork allocation and vertical align¬ 
ment optimization in road design. Hare et al. m and Lima et al. uni developed 
two mixed integer linear programming models for earthwork operation in road 
construction. Unlike Burdett et al. [3 who developed a model for the earthwork 
allocation problem considering earthwork as discrete 3D blocks, in the present 
article we use a section-based model, noting that section-based models achieve 
similar precision as 3D block based model when section lengths are less than 
30m 13 . 

While vertical alignment optimization is more complicated, it also has a rich 
research literature. In 2009 Moreb [33 developed a linear programming model 
combining the vertical alignment and earthwork allocation optimization. In 
2010, Koch and Lucet [23 advanced Moreb’s model by removing unnecessary 
errors in slope constraints. More recently, Hare et al. m incorporated the 
vertical alignment in the earthwork allocation model, resulting in a mixed integer 
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linear programming model that can be solved efficiently in practice. 

While horizontal alignment optimization is the most complicated problem, 
it has nonetheless seen a number of approaches. Jong et al. [201 ED devel¬ 
oped a horizontal alignment optimization model which was solved by a genetic 
algorithm. However, the resulting horizontal alignment does not offer any guar¬ 
antee of (local) optimality. In 2008, Easa et al. [TT] developed an optimization 
model incorporating safety constraints, which were quantified as the expected 
collisions for an alignment. Although this model guarantees global optimality, 
the associated vertical alignment cost is not incorporated in the optimization 
process. In 2009, Lee et al. [55] presented a heuristic based method to optimize 
the horizontal alignment that works in two stages. In the first stage, the heuris¬ 
tic tries to approximate a piecewise linear alignment and then in the second 
stage, it refines the solution to make the previously generated piecewise linear 
alignment compatible with a real road alignment. The solution alignment of the 
model yields a practical alignment but since a heuristic algorithm was used to 
solve the model, optimality is not guaranteed [55] ■ 

In the literature, a few studies also investigated the problem as a three 
dimensional alignment optimization problem in which the vertical and horizontal 
alignments are optimized simultaneously. 

Tat and Tao |33j proposed a three-dimensional alignment optimization model, 
which they solved using a genetic algorithm. Their model considers all of the 
major constraints in road design. Akay [5] developed a model for three dimen¬ 
sional alignment optimization for forest roads and solved it using a simulated 
annealing algorithm. Aruga |3] used a tabu search method to optimize three 
dimensional alignments of forest roads. 

A criteria-based decision support system for three dimensional alignment 
optimization was developed by Jha m considering the environmental costs. 
Jong et al. |22| presented an evolutionary model for optimizing the vertical and 
horizontal alignment simultaneously. The previous two models |I71122j were 
improved in [18j by considering accessibility, proximity, and land-use changes, 
and further improved in |25| to consider incorporating bridge and tunnel costs. 

Cheng and Lee |5] also proposed a heuristic-based model for three dimen¬ 
sional alignment optimization. The heuristic solves the models in three steps: 
first, it generates a good general horizontal alignment by adding, deleting, or 
moving the intersection points one by one, then it determines an improved hor¬ 
izontal alignment by adjusting the intersection points based on the previously 
generated horizontal alignment, and finally, it finds a better three dimensional 
alignment by tuning the vertical alignment corresponding to the previously ob¬ 
tained horizontal alignment. 

Most similar to this work, Kang et al. [24] developed a bi-level optimization 
model for road alignment design. In the upper level a set of alternative good 
alignments is generated and in the lower level an alignment is selected from 
the alternative alignments obtained in the upper level. The model [23] was 
solved using a genetic algorithm. Recently, Kang et al. [23] also proposed 
a three dimensional alignment optimization model based on genetic algorithm 
and geographic information system (CIS). 
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All of the above mentioned three-dimensional alignment optimization models 
use heuristic-based algorithms which do not guarantee optimality (or even local 
optimality), and have no or very weak convergence guarantees. 


2 Terminology 

Horizontal alignment optimization consists of finding an optimal curve connect¬ 
ing two given end-points within a designated corridor. The ground profile data 
in a corridor is given at some discrete points, named data points, within the 
specified corridor. 

There are two types of data points, namely, base data points and offset 
data points. Typically, the base data points are the points along the engineer’s 
original horizontal alignment (however, this is not strictly necessary). The offset 
data points represent the horizontal displacement from the base data points. 
The base data points are selected a few units apart between the two end-points 
along the baseline. Each of the base data points has some associated offset 
data points in both the left and the right directions, see Figure A base data 
point together with the associated offset data points is defined as a station. The 
baseline of a corridor is a curve connecting the base data points (i.e., the dotted 
curve in Figure]^. In practice, a baseline is a primarily defined alignment by 
engineers. 
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Figure 1: Corridor of a horizontal alignment. 


Base data points, offset data points, stations, and baselines are all fixed 
input data for a given horizontal alignment problem. Ultimately, the horizon¬ 
tal alignment problem is to determine the optimal curve (in terms of vertical 
alignment and earthwork cost) through those stations, subject to road design 
constraints. 

Each data point within the corridor, either a base data point or an offset data 
point, has some associated ground profile data. Therefore, for the vertical road 
profile, we can move vertically up and down for each horizontal data point. The 
horizontal and the vertical displacements from the baseline make a discrete grid 
for each station, see Figure Our goal is to find, for each station, a horizontal 
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offset that generates a horizontal alignment and a vertical alignment which is 
(locally) optimal. 



Figure 2: Three dimensional corridor. Corresponding to each horizontal offset, 
a set of vertical offsets is shown along the vertical line 


3 Geometric representation 

A horizontal alignment consists of a sequence of circular curves and tangential 
lines, which are defined by some intersection points and the radius of curvature 
associated with each intersection point. In Figure S and E are the start 
and end points of the alignment, respectively. The intersection points of the 
alignment are Pi, P 2 , and P3. Each intersection point has a radius of curvature 
that defines the circular curve. The radius of curvature associated with the 
intersection points Pi, P 2 , and P 3 are ri, r 2 , and r 3 , see Figure 



Figure 3: Geometric representation of an horizontal alignment. 
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Let i be the index of the intersection points and n be the number of in¬ 
tersection points. Since the intersection point Pi has an associated radius of 
curvature we define an intersection point with radius of curvature as {Pi, Vi), 
where Pi G and G K. Without loss of generality, we can say that the 
start and end points are two points in with zero radius of curvature and note 
them (PqiO) and (P„_i,0). So we represent a horizontal alignment PA as the 
sequence 

PA=((Po,0},(Pi,ri),(P2,r2),...,(P„-i,0)). (1) 



Figure 4: Geometric specifications of a circular curve. 

To determine the actual horizontal alignment, we need to calculate the cir¬ 
cular curves and tangential line segments from the given intersection points and 
the associated radius of curvatures. Let Pi and Fi be the left and right tan¬ 
gential point, respectively, see Figure]^ Let Ci be the center of curvature that 
corresponds to the intersection point Pi. 


3.1 Tangential point calcnlation 

We define the following variables: 6i is the angle at Pi using three consecutive 
intersection point Pi-i, Pi, and Pi+i, Qi is the intersection point of the angle 
bisector of 9i and the line joining Pi-i and Pi+i, Ui = Pi_i — Pi, V, = Pi^i—Pi, 
and Wi = Pi+i — Pi-i- We have 


9i = arccos 


U,-V, \ 

l|U.||||V,||J- 


( 2 ) 


Since the angle bisector PiQi bisects the angle 9i, we have ZPi-iPiQi = y = 
ZPi+iPiQi. The lines PiPi-i and PiPi+i are the tangent to the circle at the 
tangential point Ei and Fi. Thus we have CiEiFPiPi-i and CiFiFPiPi+i. 
The triangle APiEiCi and APiFiCi are right-angle triangles. The segment 
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PiCi is the common side of APiEiCi and APiFiCi. Since EiCi=FiCi, we have 
PiEi = PiFi. Let It be the length of PiEi. We have 


It = 


rj 

tan Y 


(3) 


Let eu- and ev^ be the two unit vector of and V^. The tangential point Et 
and Fi is calculated as 


Ei = Pi + lteu., (4) 

Pi = Pi + ( 5 ) 


3.2 Center point calculation 

Fact 3.1 ([71 Book VI Proposition III]) The angle bisector of an angle in 
a triangle divides the opposite side in the same ratio as the sides adjacent to the 
angle. 


Let If, be the length of QiPi-i. The length of PtPi-i, PiPi+i and 
are ||Ui||, ||Vi||, and ||Wi||, respectively. The line PiQi is the angle bisector of 
9i in APi_iPiPi+i. Using Fact |3.1| we have 


liu.llliw.ii 

liu.|| + |!v,||- 


( 6 ) 


Let ewi be the unit vector of W^. The point Qi can be calculated as 


Qi — Pi — 1 Ib^^^i- 


(7) 


We note = Qi — Pi- Let l^ be the length of PiCi. We have 


Let exi be the unit vector of X^. So the center point is obtained by 

Ai — Pi A Ix^'X-i ■ 


( 8 ) 

(9) 


4 Model description 

We formulate the horizontal alignment optimization problem as a bi-level op¬ 
timization problem. The inner problem is the vertical alignment optimization 
problem, while the outer problem uses this information to seek an optimal hor¬ 
izontal alignment. 
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4.1 Variables and objective function 

We can define a horizontal alignment in a specified corridor using a set of in¬ 
tersection points and the associated radius of curvatures, see Figure HI By 
varying the location of the intersection points in the associated feasible regions, 
a wide variety of horizontal alignments can be built. Note that in our model 
the feasible region of an intersection point is defined as a rectangular box. Since 
the intersection points of a horizontal alignment are in the X 2 /-plane, we denote 
Pi = {xi,yi). We write the variable vector as 

V = {{Xi,yi,ri),{x2,y2,r2) ■■■,{Xu-2,yn-2,ru-2)) 

Note that the starting point Pq and the ending point P„_i are fixed and do not 
have any corresponding radius variables. 



Figure 5: A potential horizontal alignment in a specified corridor. 

Suppose we have Us stations indexed from 0 to ris — 1. At each station, we 
define a cross-section line which passes through the leftmost offset data point and 
the rightmost offset data point. Let Uj and Vj be the leftmost and rightmost 
offset data points of a station (see Figure [^. So the cross-section lines in a 
specified corridor are written as 

Lj{t) = (1 — t)Uj + tVj for t G M, j = 0,1,..., rig — 1. (11) 

Since a horizontal alignment is a linear-circular curve, we can calculate the 
point of intersections between the cross-section lines and a horizontal alignment. 










A point of intersection can be represented as a parameter value of the cross- 
section line which is denoted as tj. If tj G [0,1] for all j, then a horizontal 
alignment is inside of the corridor. The values of tj can be mapped into the 
horizontal offsets. In our input data, horizontal offsets are discrete values, see 
Figure So we used linear interpolation to compute vertical ground profile 
data for any (continuous) offset values. The accuracy of this interpolation is 
clearly dependent on the distance between the offset data points. If higher 
accuracy is desired, then more offset data points should be included in the 
problem formulation. As computation time of this interpolation is independent 
of the number of offset data points, increasing the number of data points should 
only effect the data storage requirements, not solution time. 

Once we have vertical ground profile data at each station for a horizontal 
alignment we can solve a vertical alignment optimization problem using the 
models developed in |14j . Since we will use the vertical alignment solution, 
denoted as a black-box that we seek to minimize, it is not necessary to 

understand how exactly Cyyi is computed for this paper. Nonetheless, in Ap¬ 
pendix]^ we provide a brief overview of the mixed-integer linear program used 
to compute Cva (for curious readers). 

The final optimization problem can be written as 

min/(A) = ( ift,G[0,l], j=0,I,...,n.-l. 

' ^ ( c» otherwise, 

where Cvyi(Ar) is a function that gives the optimal vertical alignment cost for 
X. Note that when a horizontal alignment is outside of the corridor, we set the 
cost to infinity. 

4.2 Constraints 

A horizontal curve consists of tangential line segments followed by circular arcs. 
Two consecutive circular arcs are connected by a tangential line. In order for this 
tangent line to be well-defined, it is important that the circular arcs are small 
enough that they do not ‘overlap’. In Figure]^ we give an pictorial example 
of what goes wrong if the circular arcs are too large. Therefore, in order to 



Figure 6: Left: a curve with well-defined tangent lines. Right: a curve without 
well-defined tangent lines. 

maintain well-defined tangential line segments, on the line passing through the 
intersection points Pi-i and Pi, the length of Pi-iPi must be greater than or 
equal to the summation of the length of Pi-iFi-i and the length of PiEi. We 
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can write these constraints as 


+ i = l,2,...,n-2. (13) 

We will refer to these as the continuity constraints. 

If the radius of curvature is too small (or zero) then a horizontal alignment 
might get a sharp turn. So the optimal radius of curvature must be greater than 
or equal Rmin, the minimum radius of curvature. In other words, 

ri>Rmin, i = 1,2,... ,n - 2. (14) 

In our model, each intersection point has a feasible region defined by a rect¬ 
angular box, see Figure A given rectangular box is defined using the leftmost 
bottom corner point {Ixijlyi) and the rightmost top corner point {uxi,Uy^. 
Thus the box constraints corresponding to intersection point Pi are written as 

^xi ^ ^ ’^xi and ly^ Vi ^ ^ — 1 , 2 ,..., n 2 . (f^) 

5 Solution Approach 

The objective function of problem (??) is an optimization problem in itself. 
Indeed, Cva is evaluated through a large scale mixed integer linear program 
(see Appendix]^. As such, it is very hard to access the derivative information 
of the objective function (if it exists). As a result, optimizing the objective 
function cannot be accomplished by gradient based methods (such as BFGS) 
or structure based methods (such as the Simplex algorithm). Instead, it is 
necessary to apply some non-gradient based approach. We specifically turn 
our attention to two derivative-free optimization (DFO) solvers: NOMAD [5] 
(version 3.5, available at http://www.gerad.ca/nomad) and HOPSPACK [32] 
(version 2.0.2, available at http://www.sandia.gov/hopspack). In Subsections 
|5.1| and |5.2| we give very brief overviews on how each of these solvers work. For 
greater detail, we refer readers to |S| (for NOMAD) and [32] (for HOPSPACK). 

5.1 NOMAD 

NOMAD abbreviates “Nonlinear Optimization with MADS”. The solver is an 
freely available and based on the “Mesh Adaptive Direct Search” (MADS) algo¬ 
rithm designed in [1]. MADS falls into the large class of pattern search methods, 
which essentially iterate the following. 

Given an iterate a;^, a search radius > 0, and a ‘pattern’ ...d™} 

with max{|d*| : i = 1 , 2 , ...,m} = A^, evaluate f{x^ + d*) for 
each i = l,2,...m. If a better point is found, then set 

_|_ dJ and repeat using A^+^ = A^. If no better point is 
found, then set = yA*’ where 0 < 7 < 1 and repeat. 
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Pattern search methods enjoy several strong feature. First, pattern searches only 
require function evaluations, not gradient evaluations, to optimize the objective 
function. Second, unlike more heuristic methods (like genetic algorithms) in 
many cases pattern searches can be mathematically proven to converge to a 
local minimizer. Third, pattern searches are naturally parallelized. 

The MADS algorithm follows specific mathematical rules on how to generate 
the pattern, and using these rules can guarantee convergence to a local minimizer 
[4]. While the rules are complex, they are nicely described (with examples) in 

i- 

NOMAD is deterministic in nature, meaning that if the algorithm is run 
twice on the same problem, it will return the exact same outcome both times. 
As such, in our numerical tests, NOMAD is only run once per test problem. 

The NOMAD solver allows for several different stopping criterion. For our 
numerical tests, we employ the ‘minimum mesh size’ stopping criterion, which 
(in the language above) stops the algorithm when the search radius (A^) be¬ 
comes too small. 

5.2 HOPSPACK 

HOPSPACK, the Hybrid Optimization Parallel Search PACKage, implements 
the Asynchronous Parallel Pattern Search (APPS) algorithm in a C-|—I- software 
framework [271II2]. Similar to MADS, the APPS algorithm is a pattern search 
method. Although the APPS algorithm uses a different set of rules for pattern 
generation, it is also proven to converge to a local minimizer |27j . However, 
unlike MADS and APPS, the HOPSPACK algorithm includes some heuristic 
techniques to break free of local minimizers. These heuristics invalidate the 
original convergence proof of APPS, but have been found effective in practice. 
(If the heuristic are only employed a finite number of times, then the convergence 
proof of APPS still holds, but it is unclear if this safe-guard is implemented in 
HOPSPACK.) The heuristics are also stochastic in nature, so HOPSPACK will 
not necessarily return the same optimized outcome if run a second time on the 
same problem. Therefore, in our numerical tests, we have run HOPSPACK five 
times for each test problem. 

Similar to the MADS algorithm, HOPSPACK includes a stopping criterion 
based on a minimum allowable search radius (this time called ‘minimum step 
size’). For consistency with NOMAD, in our numerical tests we use this stopping 
criterion. 


6 Experimental results 

We performed numerical experiments on five different roads listed in Table 
The road profile data correspond to real roads and were given by our industry 
partner Softree Technical System Inc. 

All of the experiments were carried out in a Dell workstation with an In- 
tel(R) Xenon(R) 2.40 GHz (2 cores) processor, 24 GB of RAM and a 64-bit 
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Road Name 

Length of the Road 
(meters) 

Number of 
Intersection Points 

Road A 

358 

8 

Road B 

4000 

5 

Road C 

5013 

14 

Road D 

1156 

22 

Road E 

721 

9 


Table 1: Specifications of the test problems 


Windows 7 Enterprise operating system. Software to process and solve the 
optimization problem was implemented in C++ using Microsoft Visual Studio 
2010 Professional Edition. 

Since we are interested in a solution close to a local minimum, for both solvers 
(NOMAD and HOPSPACK) we set the stopping criterions (i.e., minimum mesh 
size and minimum step length) to 0.1. As the input data of our model are given 
in meters, the final scaling of the variables of our model goes down below 10 cm, 
which means a local minimum should exist in less than 10 cm distance (see [5] 
and [32] 1. In all tests, convergence to this level of local minimizer was confirmed 
by the algorithm. 

We used the baseline alignment of a corridor as an initial starting point for 
both solvers. Recall, the NOMAD solver gives a deterministic solution (i.e., 
different independent runs of the algorithm yield the same solution), while the 
HOPSPACK solver gives anon-deterministic solution (i.e., different independent 
runs of the algorithm might yield different solutions). So first, we solved the 
test problems using the NOMAD solver and then compare with the HOPSPACK 
solver solutions obtained by several independent runs. 

6.1 Case Study 

Before providing a summary of the results concerning all test problems, we 
focus on a single problem as a case study. Specifically, we discuss Road A, 
which is a short (358m) section of a forestry road. In Figure [Tj we visualize 
the initial baseline alignment, and the two optimized alignments (as solved by 
NOMAD and HOPSPACK). For the sake of visualization. Figure]^ also includes 
the contour lines of the terrain within the allowable corridor for Road A. 

In Figure [Tj the initial baseline alignment for Road A is presented as a solid 
black line. Notice it goes exactly down the middle of the corridor. (This is not 
a required feature of the initial alignment, it is simply because the test problem 
data, provided by Softree Inc., happened to take this form.) The optimal vertical 
alignment for the baseline alignment of Road A cost $1897. 

In Figurethe optimized horizontal alignment for Road A, when computed 
using NOMAD as the DFO solver is presented as a dashed line. Recall that 
HOSPACK in non-deterministic, and so HOPSPACK was run five times on each 
test problem. In Figure the first of these five optimized horizontal alignments 
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Figure 7: Initial and Improved Horizontal Alignments for Road A. 

Solid Black = Initial Alignment, Dashed = NOMAD optimized alignment. 
Dotted = HOPSPACK optimized alignment 


when computed using HOPSPACK is presented as a dotted line. The optimal 
vertical alignment for the horizontal alignment created using NOMAD costs 
$1361, and the optimal vertical alignment for the horizontal alignment created 
using HOPSPACK costs $1291. Both of these represent significant savings over 
the original cost. In particular, the percentage savings can be computed as 


1897- 1361 
1897 


100 % 


28.3% and 


1897- 1291 
1897 


100 % 


31.9%. 


6.2 Summary results for all test problems 

We now turn our attention to an aggregate analysis of the all test problems. 

Table shows the cost improvement of the objective functions, the number 
of black-box evaluation and wall-clock time required to solve the test problems 
using the NOMAD solver. Note that a candidate alignment that falls outside 
the corridor counts as a function call but does not require solving a MILP, which 
explains that road C and D require similar time to solve with Road D having 
double the number of function calls. The cost improvements range from 11% 
to 54% with an average of 27%. While the wall-clock times for some problems 
are high, the construction of a new road can take years. As such, the used 
computer time (which could be run over night) is not unreasonable given the 
average saving in final costs. 
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Road 

Name 

Initial 

alignment 

cost 

Optimized 

alignment 

cost 

Cost 

Improv¬ 

ement 

(%) 

No. of 
Black¬ 
box 
evalua¬ 
tions 

Wall- 

clock 

time 

(seconds) 

Road A 

1,897 

1,361 

28% 

2,073 

1,445 

Road B 

17,036 

15,198 

11 % 

2,528 

1,770 

Road C 

87,829 

69,621 

21 % 

37,165 

45,647 

Road D 

31,031 

14,418 

54% 

90,535 

47,613 

Road E 

8,054 

6,498 

19% 

11,101 

5,588 


Table 2: Cost improvement, no. of black-box evaluations and wall-clock time 
required to solve the test problems by the NOMAD solver. 


Table 1^ lists the optimum values of the objective functions and the number 
of black-box evaluations required to obtain the solution for five independent 
executions of the HOPSPACK solver for each test problem. The differences in 
the optimum objective function values are calculated with respect to the value 
obtained by the NOMAD solver. So in Table|^ a “-I-” value in the Diff. in costs 
column (i.e., difference in optimum costs obtained by the two solvers) indicates 
the HOPSPACK solver yields a better solution than the NOMAD solver and a 
“ — ” value indicates the opposite. Similarly, we also calculated the difference 
in number of black-box evaluations with respect to the number of black-box 
evaluations required by the NOMAD solver. So in Table a “ ” value the 

Dijf. in black-box eval column (i.e., differences in no. of black-box evaluations) 
indicates the HOPSPACK solver required less black-box evaluations than the 
NOMAD solver and a “ — ” value indicates the opposite. 

Combining the results obtained for the different roads listed in Table we 
make an overall comparison between the two solvers. We observed that the 
HOPSPACK solver might yield a better or a worse solution than the solution 
obtained by the NOMAD solver. Thus, considering the tolerance of the dif¬ 
ference in the optimum objective values obtained by the two solver, we count 
the number of times a solver wins with respect to the other solver. Table 
shows the comparison of the solvers for different tolerance values of the differ¬ 
ence in the optimum objective values. The x% tolerance of the difference in 
optimum objective values means if the optimum objective values obtained by 
the two solvers are in between —x% to -\-x%, then the solvers yield the same 
solution (i.e., the two solvers tie), otherwise a positive percentage value indi¬ 
cates the HOPSPACK solver wins and a negative percentage value indicates the 
NOMAD solver wins. 

In Table|^ we see that if the tolerance of difference in the optimum objective 
value is ±3% or above, then the two solver tie for more than 50% test runs (i.e., 
more than 13 test runs among 25 test runs). We can also observe that for 
any case of the tolerance change in the optimum objective function value, the 
difference in the number of times the NOMAD solver wins and the number 
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Rd 

Test 

run 

Opt. 

cost- 

NOMAD 

Opt. 

cost- 

HOPS. 

Diff. 

in 

costs 

(%) 

# 

black¬ 

box 

eval- 

NOMAD 

# 

black¬ 

box 

eval- 

HOPS. 

Diff. in 

black¬ 

box 

eval 

(%) 


#1 

1,361 

1,291 

+5.2% 

2,073 

325 

-f84.3% 


#2 

1,361 

1,423 

-4.6% 

2,073 

336 

-f83.8% 

A 

#3 

1,361 

1,418 

-4.2% 

2,073 

697 

-f66.4% 


#4 

1,361 

1,278 

-f6.1% 

2,073 

486 

-f76.6% 


#5 

1,361 

1,486 

-9.2% 

2,073 

665 

-f67.9% 


#1 

15,198 

15,510 

-2.1% 

2,528 

316 

-f87.5% 


#2 

15,198 

15,141 

-fO.4% 

2,528 

309 

-f87.8% 

B 

#3 

15,198 

15,128 

-fO.5% 

2,528 

286 

-f88.8% 


#4 

15,198 

15,172 

-fO.2% 

2,528 

392 

-f88.5% 


#5 

15,198 

15,529 

-2.2% 

2,528 

485 

-f80.8% 


#1 

69,621 

70,161 

-0.8% 

37,165 

7,547 

-f79.7% 


#2 

69,621 

70,378 

-1.1% 

37,165 

31,418 

-fl5.5% 

C 

#3 

69,621 

69,995 

-0.5% 

37,165 

3,392 

-f90.9% 


#4 

69,621 

67,301 

-f3.3% 

37,165 

3,213 

-f91.4% 


#5 

69,621 

67,045 

-f3.7% 

37,165 

4,661 

-f87.5% 


#1 

14,418 

13,190 

-f8.5% 

90,535 

15,852 

-f82.5% 


#2 

14,418 

15,154 

-7.6% 

90,535 

19,997 

-f77.9% 

D 

#3 

14,418 

14,155 

-fl.8% 

90,535 

21,194 

-f76.6% 


#4 

14,418 

14,016 

-f2.8% 

90,535 

17,816 

-f80.3% 


#5 

14,418 

15,384 

-6.7% 

90,535 

21,339 

-f76.4% 


#1 

6,497 

6,524 

-0.4% 

11,101 

2,019 

-f81.8% 


#2 

6,497 

6,475 

-fO.3% 

11,101 

4,332 

-f60.1% 

E 

#3 

6,497 

6,497 

0 .0% 

11,101 

1,996 

-f82.0% 


#4 

6,497 

6,502 

-0.1% 

11,101 

2,501 

-f77.5% 


#5 

6,497 

6,476 

-fO.3% 

11,101 

3,120 

-f71.9% 


Table 3: Comparison of optimum objective function values and required no. of 
black-box evaluations for the HOPSPACK and NOMAD solvers to solve the 
test problems. 


of times the HOPSPACK solver wins is at most 2. So in terms of optimum 
objective values obtained by the two solvers, the performance of both solvers 
are roughly equivalent. In Table we can see that for all of the 25 test runs, 
the HOPSPACK solver required on average 78% less black-box evaluation than 
the NOMAD solver. So the HOPSPACK solver is roughly five time faster than 
the NOMAD solver to compute the optimum solution. 
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Tolerance of the 
difference in the 
optimum costs 

No. of times 
the NOMAD 
solver wins 

No. of times 
the 

HOPSPACK 
solver wins 

No. of times 
the two 
solvers ties 

± 1 % 

8 

7 

10 

± 2 % 

7 

6 

12 

± 3% 

5 

5 

15 

± 4% 

5 

3 

17 

± 5% 

3 

3 

19 

± 6 % 

3 

2 

20 

± 7% 

2 

1 

22 

± 8 % 

1 

1 

23 

± 9% 

1 

0 

24 

± 10 % 

0 

0 

25 


Table 4: Overall comparison of the HOPSPACK solver and the NOMAD solver 
with the optimum objective function values. 


7 Conclusion 

Determining a good horizontal alignment is the first concern when designing a 
new road. In this paper, we model horizontal alignment using the geometric 
specifications which are used by engineers in practice. We further pursued 
a novel approach to address the horizontal alignment optimization problem. 
While most of the studies in the literature used heuristic based methods, we used 
a derivative-free optimization approach since it converges to a locally optimum 
solution. Our method obtains a horizontal alignment that is locally optimal and 
a vertical alignment that is globally optimal. On our road samples, this locally 
optimal solution was on average 27% cheaper than the solution provided by the 
engineers. As the solution of the optimization model yields a practical horizontal 
alignment that satisfies geometric specihcations and engineering requirements, 
the results represents signihcant savings in the final design of a new road. 

Backward bends in a horizontal alignment (i.e., roads where the cheapest 
path from start to end requires moving away from the target destination for some 
portion of the horizontal alignment) are known to be a challenge to optimize 
[301 page 123, Chapter 5], m page 21, Section 2.4.3], |3T]. (Such behaviour 
is often referred to as backtracking.) However, our model can generate such 
alignments without difficult. 

Although our proposed model works well for solving practical horizontal 
alignment optimization problems, it can be improved further for better preci¬ 
sion and performance. In our model formulation, we considered that the cross 
sections in a corridor are hxed and taken corresponding to the baseline align¬ 
ment. However, cross sections should be always perpendicular to a horizontal 
alignment. When a horizontal alignment is significantly different from the base¬ 
line alignment, a set of new cross sections should be generated to increase the 
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precision before calculating the corresponding vertical road profile. 

A surrogate cost function is an approximation of the original cost function 
which is cheaper to compute. Using a surrogate cost function might reduce the 
solution time required to get a solution. (The NOMAD solver can exploit the 
usage of a non-adaptive surrogate.) 

In our model, we only consider the construction costs to formulate our cost 
function. In the future, land acquisition costs could be incorporated by con¬ 
sidering the unit cost of a piece of land corresponding to two consecutive cross 
sections in a corridor. We can also include pavement costs. 

Both solvers need an initial starting point to start the algorithms. Since 
solutions obtained by both solvers are locally optimum, we can use multiple 
starting points (i.e.; multiple initial alignments) to obtain a better solution 
quickly. How to choose multiple alternative good alignments in a specified 
corridor is the subject of active investigation. 

During the optimization process, at each iteration, both of the derivative 
free optimization solvers solve a large number of vertical alignment optimiza¬ 
tion problems (i.e., large scale mixed integer linear programming (MILP) prob¬ 
lems). At the earlier stage of the optimization process (i.e., when the mesh size 
is coarse) we can relax some of the parameters of the vertical alignment opti¬ 
mization problem to get an approximate cost and then at the later stage (i.e., 
when the mesh size becomes relatively small) we can tighten the parameters 
to get accurate costs. This policy might reduce the solution time signihcantly. 
We can also use a warm start of the vertical alignment optimization problem 
when the horizontal alignments are close to each other to accelerate the vertical 
alignment optimization process. So the interconnection between the derivative 
free optimization solver and the MILP solver can be a potential way to reduce 
the solution time. 
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A Vertical Alignment Optimization Model 

In this appendix we give a very brief overview of the Quasi Network Flow (QNF) 
MILP used to solve the vertical alignment optimization subproblem. The model 
is fully described in [M] and we refer there for details and definitions. The 
following equations are only provided to convince the reader that the model is 
indeed a MILP. 

The variables are , 1^^“, ui, and binary variables yk,t- The 

objective function is 


min PzV++ Y + 

iG5uB iG5UW iGS 

teT 

J6B jew 

ter ter 


The model includes numerous linear constraints like conservation of material 


rr j_ fu 


E 


j,i+l,t 


leB 

■d{j)=i 


rr _[_ ru 

' /i,z—l,t 
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f° 

■t j,i--i,t 


j&B 

^(j)=i 




E 




iew 


for alH G T, * S S, 


fli-i,t + fi+i,i,t 


+ E /r+i 


■.j,t ! 


jew 


for all t G T,i G S, 


and volume balance constraints 

Yf"^-pt+m+i,t+m,t = v^, 

ter 

E fUl,i,t + fi,i,t = > 


ter 


j,‘s{j)-^,t fj,‘&{j)+l,t ^j ) 

ter 

E + fipij}+i,j,t — > 

ter 

The volumes are related to the quadratic spline 


Vi G 5, 
Vi G 5, 
Vj G B, 
Vj G W. 


P{s) 


'A(s) 

P2is) 


if '^5(2,1) < 5 < 5^(2, 


if '^5(5,1) — ^ — "^( 5 ( 1 ,rig)' 
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where Pg{s) = Ug^i + ag^2S + o,g^^s^ by the offset value 


The smoothness of the quadratic spline is enforced by 

Pg-l{ss{g,l)) = Pg{ss{g,l))t '^9 ^ G \ {!}) 

Pg-As5ig,l)) = PgisSig,!)) Vffeg\{l}, 

while grade constraints involves the derivative of the spline 

GL<P'issig,i))<Gu, V5Ge\{l}. 

Other constraints include hxed point constraints, side-slopes constraints, 
block constraints, and bound constraints. All these constraints are linear but 
may involve binary variables. We refer to M for their explicit definition and 
expression. 
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